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Abstract The Cassini spacecraft collects high resolution images of the saturnian satel- 
lites and reveals the surface of these new worlds. The shape and rotation of the satellites 
can be determined from the Cassini Imaging Science Subsystem data, employing limb 



coordinates and stereogrammetric control points. This is the case for Epimetheus (Tis- 



careno et al.||2009 l that opens elaboration of new rotational models (Tiscareno et al 



2009 Noyelles[ 2010 Ro butel et al.|20lT| . Especially, Epimetheus is characterized by its 



horseshoe shape orbit and the presence of the swap has to be introduce explicitly into 
rotational models. During its journey in the saturnian system, Cassini spacecraft accu- 
mulates the observational data of the other satellites and it will be possible to determine 
the rotational parameters of several of them. To prepare these future observations, we 
built rotational models of the coorbital (also called Trojan) satellites Telesto, Calypso, 
Helene, and Polydeuces, in addition to Janus and Epimetheus. Indeed, Telesto and 
Calypso orbit around the L4 and L5 Lagrange points of Saturn- Tethys while Helene 
and Polydeuces are coorbital of Dione. The goal of this study is to understand how 
the departure from the Keplerian motion induced by the perturbations of the coorbital 
body, influences the rotation of these satellites. To this aim, we introduce explicitl y the 



perturbation in the rotational equations by using the formalism developed by (Erdi 



19771 to represent the coorbital motions, and so we describe the rotational motion of 



the coorbitals, Janus and Epimetheus included, in compact form. 
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1 Introduction 



The space mission Cassini orbiting Saturn since 2004 provides a lot of data and notably 
images of the surface of over 20 satellites (Thomas 2 010[). From th ese images methods 
based on limb profile provide the satellites shape ( jThomas "et~aL|T998l [2007) |Thomas| 



20101. To combine the different images taken at different time, it is crucial to assume 
a rotational motion of the satellites. Thus, the rotational determination is obtained 
as a by-product result. Usually for almost satellites a synchronous uniform rotational 
motion seems to be enough to fit the images | Thomas et al.|2007 1. However, for specific 
cases such as Janus and Epimetheus more elaborated rotational models are required 
( Tiscareno et al.|2009 1 because of the swap in their orbital motion. This swap is due to 
the 1:1 mean motion resonance where every 4 years these satellites approach and swap 
their orbits by a few tens of kilometers. The resulting orbits in the rotating reference 
frame are an horseshoe shape orbits. 

In the Saturnian system, four additional coorbital satellites (i.e. in 1:1 orbital reso- 
nance) are currently known. They are Calypso and Telesto that are coorbital of Tethys, 
and, Helene and Polydeuces, coorbital of Dione. Contrary to Janus and Epimetheus, 
their orbits describe a tadpole shape in the rotating reference frame (Fig. [T]). Their 
masses are very small with respect to their main satellites Tethys and Dione, whereas 
Janus and Epimetheus have a comparable size. According to Horizons ephemeris 
(Gio rgini et al.|["l996 l, the mass ratio between Janus and Epimetheus is about 3.604, 
while the masses of Telesto, Calypso, and Helene are several tens of thousand times 
smaller than the one of Tethys and Dione. We note that Polydeuces' mass is currently 
unknown, but it is expected much smaller than the one of Helene because the mean 
radius of Polydeuces is less than 2 kilometers (Porco et al. 20071. Consequently, the 
orbital motion of these four small satellites may be deduced from the restricted three 
body problem including Saturn and one of these small bodies, Tethys or Dione. In this 
framework, the small satellites are located close to the L4 and L5 Lagrange points 
of the main satellites, executing a tadpole orbit around one of the Lagrange points 
(see |Christou et al.| ( |2007| ) and references therein). The objective of the present paper 
is to investigate the influence of the orbital oscillations resulting from the 1:1 orbital 
resonance, on the rotational motion of these satellites by developing a general method 
based on the Hamiltonian approach and parametrization of the orbital motion. 

The coorbital satellites have irregular shape that probably results from their orig- 
inal accretion or their impact history (Porco et al. 2007 Thomas 20101. The shape 
can be approximated by an ellipsoid obtained from the best fit to the numerical shape 
model developed by Thomas (2010). But, for the rotation, the main parameters of 
interest are the moments of inertia. They can be deduced from the semi-axes of the 
ellipsoid by assuming an homogenous interior. The error can be estimated of the order 
of 30 % for Janus (see Section [4j| , so we could explore a large range of values around 
the homogenous ellipsoidal shape. In addition, we assume that the rotation of these 
satellites is synchronous as in |Porco et al. (20071 and 



Thomas 



(2010) because this is 



the most expected state due to tidal dissipation (e.g. Peale 19991. 

In Robutel et al. (2011 1 a perturbative approach has been used to solve the equa- 
tions of the rotational motion by assuming that the orbital trajectories of Janus and 
Epimetheus were quasiperiodic. This method allows to determine and to measure the 
influence of the coorbital resonance on the rotation of theses satellites. Here, we propose 
to apply such a method to compute efficiently the rotational motion of these satellites. 
As the satellites present a 1:1 orbital resonance and 1:1 spin-orbit resonance, first we 
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clarify vocabulary used in this paper. On one hand, the 1:1 mean motion resonance 
generates tadpole or horseshoe orbits along which the difference of the longitude of 
the two satellites are perturbed (called usually libration) with an orbital frequency 
v that we called the coorbital frequency. On the other hand, the 1:1 spin-orbit reso- 
nance also causes a libration corresponding to an oscillation of the rotation angle of the 
satellite with respect to an uniform rotation. For quasiperiodic motion, the frequency 
that depends on the amplitude of the oscillation is called proper libration frequency. 
This frequency will be denoted a at the center of libration (see e.g. Wisdom ( 2004 1). 
The term libration is used only to describe the quasiperiodic rotational motion of the 
satellites. 

The paper is divided as follows: first we describe the representation of the orbital 
motion that we use to model the coorbital satellites. We then present a Hamiltonian 
formulation of the rotational problem and discuss on non-resonant and secondary res- 
onant cases. These secondary resonant cases result in commensurability between the 
libration frequency and orbital frequency. They appear for some values of the triaxiality 
of the satellites and we detail the 1:1 and 1:2 secondary resonances. In particular, we 
show how the overlap of resonances associated with the coorbital frequency generates 
relatively large chaotic regions. Finally, we discuss our results in prevision of future 
observations and the assumption made in this paper. 



2 Approximation of the orbital motion 



In this paper, we restrict the orbital model to the circular restricted three-body prob- 
lem. It is within this simplified framework that the phenomena we want to highlight 
will appear the most clearly. Let us consider two bodies of masses mo (Saturn) and 
mi (Dione or Tethys) with mo 3> m\ and 



M = 



m\ 



mo + mi 



(1) 



The body of mass mi describes a circular orbit of radius ro centered on the most 
massive body with an angular velocity (mean motion) equal to n. The origin of the 
reference frame is located at the massive body (the planet), and polar coordinates are 
ass igned to the two other bodies: (ro,v p ) to the main satellite, and (r,v) to the trojan. 



Erdi 



( 19771, in his paper dedicated to the dynamics of the Jovian Trojans, derives 



asymptotical solutions of the trojan motion in the case of the elliptic restricted three 
body problem expending the series up to order two in -JJi. Limiting Erdi's theory to 
the first order in y^Z and assuming circular revolutions for the two massive bodies, the 
coordinates (r, v) can be approximated by these expression^] 



r(t) 



ro 



1 + y/Ji pcos(v p (t) + C(t) + ip) 



2 d£(t) 



3n dt 

v(t) = v p (t) + C(t)-2y/J2psm(v p (t) + C(t) + iP), 
v p (t) = nt + Vpfi, 

df 2 ' 



i-(i) = -3rm 2 (l - (2 - 2cosCWr 3/2 ) sinC(t). 



(2) 

(3) 
(4) 

(5) 



1 See the appendix for the relationship of these expressions with the Erdi s theory 1 1977 I 
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The coefficients (p, ip, v p fi) are constant parameters linked with the initial conditions. 
This model provides a quasiperiodic approximation of the motion of the massless body 
in the coorbital 1:1 mean motion resonance, which is valid for a tadpole orbit as well 
as for a horseshoe obit. 

The previous expressions split in two components: a fast component whose ampli- 
tude is controlled by p, and a slow component governed by the function This function 
represents the temporal variations of v — v p averaged on one orbital period 2n/n. An 
usual approximation of this quantity is given by the solution of the differential equa- 
tion ([5} (see |Kevorkian| ( |1970| ), |Salo and Yoder| ( |1988[ >, and also |Morais| ( |2001| ) for an 
averaged Hamiltonian formulation). This equation provides a good approximation of 
the coorbital secular motion as long as the trajectories stay outside of the Hill sphere 
(here a circle of radius ro(p/3) 1 ^ 3 whose center is the body of mass mi). According to 



Kevorkian ( 1970 1, this condition imposes the parameter p not to be very large with re- 



spect to unity. The system associated with the differential equation (|5| possesses three 
fixed points corresponding to the two stable equilibria L4 and L5 for (£, £) = (±7r/3, 0), 
and to the unstable equilibrium L3 when (C, C) — (""jO). The collinear points L\ and 
1/2 ar e not fixed points of the equation (jSJ, but as the points lie on the boundary of 
the Hill sphere, we do not consider this problem. In contrast to the real problem, in 
this one degree of freedom approximation given by the equation |5|, the stable and 
unstable manifolds associated with L3 coincide (the system is integrable). The domain 
inside these manifolds is filled with tadpole orbits: periodic orbits surrounding L4 or 
-L5 whose frequency varies from y / 277m/2 in the neighborhood of the corresponding 
equilibrium to approaching the separatrix. Outside of the separatrices lie the horse- 
shoe orbits whose frequency ranges from zero on the separatrix to infinity when £ tends 
to infinity. However, the value of the frequency is bounded because when [f| is large, 
the distance between the trajectory and the collision with the mass mi is small, that is 
contrary to our hypothesis that the horseshoe trajectory must not enter the Hill sphere 
(see Eq[5]|. Therefore, it can be shown that the frequency of the horseshoe trajectories 
is bounded by a quantity of the order of y/pn, similar to the frequency in the tadpole 
regime. 

The fast component of the motion is equal to ro^/Jlp cos(u p + f(t) + 1J1) in the 
expression of r and equal to — 2^/~ppsm(vp + £(i) + ip) for v. These two functions 
are quasiperiodic with frequencies n and v where v — O(^fp) is the frequency of the 
periodic function £ representing the secular motion. If we do not consider the function 
C, the expressions of r and v are the same as those obtained by a first order expansion 
in eccentricity in the Keplerian case, the eccentricity being equal to ^fpp. This remark 
will be exploited below. 

For the real trojan satellites of Saturn, formulas |2]| to (|5| provide a quite rough 
approximation of their orbital motion. Indeed, even if we can assume that the motion 
of Dione (resp. Tethys) is not too much affected by their coorbital companions Helene 
and Polydeuces (resp. Telesto and Calypso), the motion of Dione and Tethys is not Ke- 
plerian. In particular, the 2 : 1 mean motion resonances between Dione and Enceladus 
and between Mimas and Tethys generate long period variations of the orbital elements 



of these bodies (period of about 70.5 years for Tethys and 11.1 years for Dione (Vienne 
and Duriez 19951). Here, the effect of Saturn's oblateness is not directly taken into 



account. Its main impact on the orbital motion is to slightly modify the fundamental 
frequencies. As we derive these frequencies from the Horizons ephemerides that already 
include several spherical harmonics related to the shape of Saturn, the oblateness is im- 
plicitly considered. In addition, a second effect of the Saturn J2 factor is to modify the 



5 



triangular equilibrium configuration. Indeed, the equilateral triangle flattens slightly 
to give an isosceles triangle (see Sharma and Rao 19761, but this modification is very 
small in regard to the other neglected perturbations. Even if the orbital model that 
we consider here is not accurate enough to describe precisely the orbital motion of the 
studied satellites, its injection in the equations of the rotation is enough to describe the 
rotational libration with a better accuracy than the accuracy providing by the Cassini 
data ( jTiscareno et al.|2009 l. 

The parameters used in this paper to model the orbits are displayed in the table[T] p, 
is the mass of the main satellite associated with the coorbital in the three-body problem 
divided by the mass of Saturn. The value of p is fitted to the mean eccentricity of the 
satellite, derived from the ephemerides, because it controls the amplitude of the short 
time variations of r and v — nt therefore it acts like the mean eccentricity of the satellite. 
This crude approximation is sufficient to give a reasonable order of magnitude of the 
libration amplitude for the coorbital satellites (see Table[2|. The secular variations of v 
and r are represented by the function £, solution of the differential equation |BJ. For the 
first four satellites of the table, this solution is chosen such a way that its period 2n/v 
(y is given in the last column of the table) is the same as the period of the libration 
around L4 or L5 (depending on the satellite) deduced from the Horizons ephemerides 
( jGiorgini et~aL][ll96| ). 

For Janus and Epimetheus, the situation is more complicated because the long- 
term component of their true anomaly (the average of v with respect to the orbital 
period, namely the 8-years component) does not verify the equation djjl. But, as it is 
shown in Robutel et al. (20111, the average of the relative mean anomaly (difference 
between Janus' anomaly and the Epimetheus one) satisfies the differential equation (JSJ. 
In order to describe the long term temporal evolution of r and v we use the expressions 
@ and j3j again, but by replacing the function £ by m^;/(mj + rriE)C f° r Janus, 
and mj/(mj + m^C, for Epimetheus, where ( represents the above mentioned long 
time variations of relative orbits of Janus and Epimetheus, is solution of the equation 
|5|. More explanations concerning this point can be found in Robutel et al. (20111 



and references therein. The minimum and maximum values of £ are presented in the 
fourth and fifth columns of Tab. [T] For the orbital motion of Janus and Epimetheus, 
the variation in ( is greater than 180 degrees resulting from the horseshoe shape orbits 
of these satellites. On the contrary, the variations in £ are smaller for the satellites in 
tadpole orbits. 



satellite p p min£ mm( n v 



(deg.) (deg.) (rad/day) (10 3 rad/day) 



Polydeuces 


2.910" 


-6 


11.5 


268 


321 


2.3 


7.93 


Helene 


2.910" 


-6 


1.(5 


47 


77 


2.3 


8.18 


Telesto 


1.010" 


-6 


0.6 


296 


304 


3.3 


9.03 


Calypso 


1.010- 


-6 


0.8 


59 


61 


3.3 


9.02 


Janus 


9.710" 


10 


218 


6 


354 


9.03 


2.15 


Epimetheus 


3.310" 


-9 


171 


6 


354 


9.03 


2.15 



Table 1 Orbital parameters of the six coorbital satellites deduced from the Horizons 
ephemerides: p = m\/(m\ + mo), p that appears in the expressions |2| and |3| is derived 
from the satellite's averaged eccentricity, n is the mean motion and v the coorbital frequency. 
The fourth and fifth columns give minimal and maximal values of the long-term component f 
satisfying the equation d6l. 
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Fig. 1 Average orbit of the satellites modeled by the equation 15k plotted in the plan 
(£,d£/dr). The coordinates of L3, L4 and L5 are respectively (tt, 0)^(n/3, 0) and (57r/3,0). 
The two tadpole orbits surrounding L4 are associated with Calypso (smallest one) and Helene, 
while the ones around L5 corresponds to Telesto (smallest one) and Polydeuces. The outermost 
orbit describes the average motion of Epimctheus with respect to Janus. 



In order to display the average orbits of the 6 coorbital satellites on the same plot, 
different time-scale is used for each body. Rather than t, we use the time r = y/Jlnt in 
such a way that the differential equation Q becomes free from any parameter, that is 
to say: 

= -3(l-(2-2cosCr 3/2 )sinC (6) 

This rescaling allows us to compare these orbits in the plan (£,d(/dr) and to add 
the separatrix associated with L3 (black curve) which discriminates the tadpole orbits 
surrounding L4 and L5 from the horseshoe orbits outside of the separatrix. For the 
four satellites in tadpole orbit (Polydeuces, Helene, Calypso, Telesto) the rescaling 
coefficient y^n is derived from the table [TJ while the initial conditions of at r = 
are [C,max, 0) where C,max is given in the fifth column of the table[l] As mentioned above, 
the individual orbit of Janus or Epimetheus is not directly related to the equation dHl, 
so we use their relative orbit (red in figure [T|. This orbit is the integral curve of |6| 
with r = \/ (mj + m^)/ {mg + mj + m^)nt which includes the point of coordinates 

{Cmax, 0). 
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3 A Hamiltonian formulation of the rotation 

3.1 A simple quasiperiodic model 

The dynamical equation governing the physical libration is inferred from the angular 
momentum balance equation between the body's inertia and the gravitational torque 
exerted by Saturn projected onto the equatorial plane of the body. These equations 



have the following form (see e.g. Robutel et al. (2011 1) 



, sm2(6-v) = 0, wither z = v - ' , (7) 

r J rg C 

where 6 is the rotational angle, G the gravitational constant, tuq the mass of Saturn, 
A, B, C are the principal moments of inertia such as A < B < C, and (r, rrj, v) are the 
orbital parameters defined in the previous section. The frequency a is constant and 
corresponds to the frequency of libration at the equilibrium point, a depends on the 
difference B — A and for a spherical body, or an axisymetric body along the north-south 
direction, it is equal to zero. By introducing the angular distance x = 9 — v between 
the long axis of the satellite and the direction satellite-planet, we obtain the equation 



(—} 3 sin 2x = -2^Jj,pn 2 sin(v p + £(t) + tp), (9) 



* + T \ r) Sin2X = (8) 

whose solutions describe the motion of the planet in the satellite's sky. 
Using the approximation of v given by the relation {3b, we get 

f ' ' ' 

where the terms of order fj, are neglected. As mentioned in section [2] r and the right 
hand side of the equation (9| are quasiperiodic functions of frequencies n and v, thus 
the differential equation (|9l depends quasiperiodically on the time. 

The mass ratio p being small, we adopt a perturbative method. We begin by as- 
suming that fi — 0. The equation which is now equivalent to the one representing 
the motion of the pendulum, possesses a stable fixed point for x = x = surrounded 
by periodic trajectories whose frequencies are close to a in a small neighborhood of 
the fixed point. Now, increase the value of \i corresponds to perturb the pendulum 
equation in a quasiperiodic manner. If the size of the perturbation (governed by /x) is 
small enough, and under other suitable conditions excluding some resonances between 
the frequencies cr, n and v, the fixed point of the unperturbed equation is replaced 
by a quasiperiodic solution of frequencies (n, v) around which revolves quasiperiodic 
orbits possessing 3 frequencies (see |Jorba and Sim6| ( |1996[ )). The frequencies of the 
central orbit are the two forcing frequencies n and u, while the other orbits possess an 
additional "free" or "proper" frequency (close to a) which depends on their initial con- 
ditions. The perturbed systems thus possesses a one parameter family of 2-frequencies 
quasiperiodic orbits which end in x — x — when p tends to zero. 

These forced solutions play a major role in the rotational dynamics of the close 
satellites. Indeed, in presence of dissipative phenomena, the forced solution will be a 
stable equilibrium (quasiperiodic attractor) towards which most of the trajectories will 
converge (see |Celletti and"Ch icrchia (2008)). It is worth mentioning that the quasiperi- 
odic attractor of the dissipative case and the quasiperiodic forced solution of the con- 



servative case are not the same trajectories (see Landau and Lifchitz 1960 chap. 28) 



s 



But their frequencies, which are the orbital frequencies, are independent from the dis- 
sipation and consequently are equal. Let us also mention that for small dissipations, 
the two trajectories are very close, and that the dissipative limit cycle tends to the 
conservative forced orbit when the dissipation vanishes. For these reasons, we will de- 
vote the rest of this article to the study of the dynamics in the neighborhood of the 
forced quasiperiodic orbit. 



3.2 Hamiltonian formulation of the problem 
3.2.1 Autonomous Hamiltonian system 

In order to adopt an autonomous Hamiltonian formulation which takes into account 
the different time-scales, let us denote: Ai = v p + ip and A2 = y/Jint. The solution £ of 
the equation ([5| can be rewritten as: 

CW = C(A 2 ^ 1 n- 1 ) = C(A 2 ) (10) 

and its time derivative as: 

^(t) = ^nC'(A 2 ) J with c'(A 2 ) = ^(A 2 ). (11) 
It turns out that: 

=l + vWi,A 2 ) + 0(|i), where 
h(Xi, A 2 ) = 2£ (A 2 ) — 3/9COSW, with u = Ai + £(A 2 ) 



and 



i)(t) =ra(l + Vw(Ai,A 2 )) + 0(R) I where 
ffCAi, A 2 ) = C'(A 2 ) - 2pcosit. 



(13) 



Finally, using the conjugated variables (x, X) = (6— v(t), 9 — n), (Ai, Ax) and (A 2 , A2), 
the equation |9f is equivalent to the canonical equations associated with the Hamilto- 
nian H given by: 

X 2 a 2 

H = nAi + ^J2nA 2 + — - — (1 + <JfIh(\i, A 2 )) cos2x - n^Jig{\i, A 2 )X (14) 

As we are only looking for the solutions in a neighborhood of zero, we expand the term 
cos 2a; at order 4. Using the variables (y, Y) = (a 1 ^ 2 x,a^ 1 ^ 2 X), the expansion of the 
previous Hamiltonian (that we still denote by H) reads: 

y 2 + Y 2 
H=nA 1 + (T— ^ 

/ / 2 4\ \ 0-") 

+ V^(nA 2 + f o V — - V — \ ^(Ai,A 2 )-n V / o : ff(Ai,A 2 )y 

In this expression the terms depending only on Ai and A 2 are omitted. 

An usual way to study the dynamics for small values of y and Y is to reduce the 
Hamiltonian (151 to a normal form (see e.g. Arnold et al. (2006l). This will be done 
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in two consecutive steps. The first one, performed in Section [3.2.2| leads to get rid of 
the linear terms in Y, that is to eliminate g(\i, X 2 )Y using a new set of canonical 
coordinates. This is equivalent to bring back the force trajectory to the origin of the 
Cartesian coordinates. The second step, which is described in Section |3.3| consists in 
the replacement of the Cartesian coordinates by polar symplectic ones before averaging 
the obtained expression. 



3.2.2 The quasipenodtc forced trajectory 



The suppression of the term g(Xi, X 2 )Y in the expression ( 15 I is performed by in- 
troducing a new coordinate system (v , Y , X'i, A±, X' 2 , A 2 ) linked to the coordinates 
(y, Y, Ai, A%, A2, ^2); by the relational] 

, OF v i^dF , dF A ' ^dF 

y = y -gY^ Y = Y + w Xj = Xj ~m Aj = Aj + W (16) 

where the generating function F is a solution of the partial differential equation: 



DF ( dF dF \ — 
l dX[ + a y dy 7 ~ V 'd~Y') = n ^Vw?( A 'l> \' 2 )Y'. (17) 



We choose here the solution: 

F = f^y' - 2p ^cos(A- 1+ C(A 2 ))-^sin(A- 1+ C(4)). ( , ^ ( 



This solution is valid only if a 7^ n, and more generally when \a — n\ is not too small, in 
other words when the system is far enough from the 1:1 secondary resonance between 
the libration frequency a and the mean motion n. 

In the coordinate system (y' , Y' , X'i , A\ , X 2 , A' 2 ) the Hamiltonian (that we denote 
H again) becomes: 

12 1 yl2 li / / 12 /4 \ \ 

H = nA[ + a V Y + (nA' 2 + (o- V — - V —\ h(X[, A 2 )J . (19) 

It is now clear that the trajectory 

(y'(t), Y'(t), Ai(*), A 2 (*)) = (0, 0, nt + X'^^nt + X' { 2 ] ) (20) 



is a quasiperiodic solution of (19 1. This trajectory depending only on the orbital fre- 
quencies is the forced orbit which reads in original coordinates: 

,2 



x(t) = 6(t) - v(t) = -2p^^^ sin(A' 1 + C(A 2 )). (21) 



The expression of the forced trajectory in a reference frame rotating with the angular 
velocity n reads: 

9{t) - v p (t) = CM - 2 Pv ^ f 2 sin(Ai + C(A' 2 )). (22) 



2 Let us notice that the transformation Q 16p is canonical only when the terms of order two 
and more in ^fji are neglected 
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The solution is the sum of two terms of different nature: a long-period term, which 
amplitude is about 1 



for Calypso and can reach 135° for Epimetheus (Robutel et al 
20111, and a term at the orbital period whose amplitude, given in Table |2j does not 



exceed 9 . The amplitude of this libration is very similar to the one obtained for a 
satellite in Keplerian motion. If we equate the term p^/Ji with the orbital eccentricity, 
we recover the "Keplerian" amplitude expression but the phase has a very different 
behavior. In contrast with the Keplerian case, the phase is a quasiperiodic function 
of the time and the solution is close to the solution obtained by adiabatic invariant 
method developed in Robutel et al. (20111. 



3.3 The non-resonant case 



In order to proceed to the reduction of the Hamiltonian (19 1, let us introduce the 
symplectic polar coordinates 

(y',Y') = (V2Wsmw,V2Wcosw). (23) 



In these variables, the Hamiltonian ( 19 I can be split in three components: the unper- 
turbed part Hq given by 

H = nAi + aW, (24) 
the averaged part of the perturbation Hi which reads 

— / W 2 ( W 2 \ ~, , 

H 1 =n^A' 2 - — +^laW- — U'(\' 2 ), (25) 



and its remainder Hi given by 



W 2 \ „ „ _ ;/ , W 2 



Hi = {-^aWf + (1 + 2^i(')—j Qu - (1 + 2^C')^<?o,4 



(26) 



with 

2C p .q = cos(ptt + qw) + cos(pu — qw) — 2 cospu cos qw (27) 

where u = + CC^a)- 

These expressions show that the angles w and Aj vary rapidly (ib = 0(a) and 
X'l — 0(n)) while \' 2 is a slow angle (A2 = O(^fjin)). Consequently, an usual way to 
find the approximate solutions of the previous Hamiltonian system is to consider its 
average with respect to the fast angles Ai and w. 

If, as it is the case in the previous section, the terms of order p, and more are 
neglected, and if the term cos 2x is expanded in Taylor series up to degree 2p in x, the 



resonances that occur during the first order averaging process of the Hamiltonian ( 26 1 
are n + ka = 0, k being an integer satisfying \k\ < 2p. 

As mentioned above, the Hamiltonian H can be reduced to its average with respect 
to the fast angles A^ and w given by the expression: 

/ , W 2 ( W 2 \ -1 1 

K = nA'i +aW + n^A' 2 - — + s/Jl laW - — J C (A' 2 ), (28) 
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where the terms which are negligible with respect to ^fji and W 2 (when ji and W 
tend to zero) are omitted. The transformation linking H to K is approximated by the 
expression: 

dF a „ dF a , , 

z = z ~-dz' z = z+ ^¥' (29) 

where (z, Z) is a couple of conjugated coordinates belonging to {(A' l5 yl^ ) , (\' 2 , A 2 ), (w, W)}, 
(z, Z) are the average variables, and F a the generating function of the transformation. 
For simplicity, the "bar" symbols will be omitted in the rest of this section. 
The function F a being a solution of the equation 

dF a dF a 2- , . 

n^i-+a— \-Hx = 0, (30) 



we adopt here: 



with 



d\i dw 



Fa = -Y, h pAW, >n)S p , q , (31) 

P,Q 



if am(pu + qw) sm(pu -qw)\ 

2 \ pn + qa pn - qa J ^ v (32) 

Sp t q — if H Pt q = 0. 

Solving the Hamiltonian system associated with K allows to give the solutions of 
the problem in the original variables [x, X). In particular, we have: 



,2 



x(t) — ~1JJi,p—Fp^ — =■ sin u + \/ 2Wa ^ sin w 

" U 34 4 ^ 

+ V h w [? (A2) sinTt + S (SS sin( ^ +w) + ht sin( ^ - *>) ] ' 

where is a constant and the temporal expressions of the angles u and w are: 

u(t) = nt + C(*)+^, 

( W \ ( 34 ) 
to(t) =w + I a - — 1 t+ C(t). 



When W — 0, we find the forced trajectory (21 1 whose fundamental frequencies are 



only the orbital frequencies (n, v). When W > 0, the trajectories which oscillate around 
the forced trajectory have one additional frequency that depends on the distance to 
the forced orbit W: a(W) = a - W/2 + 0(W 2 ). We note that, unlike the naive 
representation which consists in separating the solution of the motion in a forced and 
a free component, the part of the trajectory which depends on the parameter W does 
not contain only a(W) but mixed term containing the three main frequencies. 
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3.4 Secondary resonances 



3.2.2 



that due to the presence of the denominator a A 



We have seen in Section 
the transformation defined by the relations (161 to (18 1, which maps a neighborhood 



of the forced orbit onto a neighborhood of y = Y = 0, is not denned when the system 
approaches the 1:1 secondary resonance. Similarly, the averaging performed in Section 
3.3| makes sense only when n + ka ^ 0. Consequently we will see on two particular 
cases (the 1:1 and 2:1 secondary resonances) how the topology of the phase space is 
modified when secondary resonances occur. Such a development of secondary spin-orbit 
resonance has been done by |Wisdoni] ( |2004j ) in the case of Keplerian orbit and more 
specifically for the 3:1 resonance. 



3.4-1 The 1:1 secondary resonance 



In order to study the dynamics in the neighborhood of 1:1 resonance between the 
frequency of libration and the orbital frequency, we start from the initial Hamilto- 
nian (151. Using the new coordinate system (w, W, A,-, A A defined by the relation 



(y, Y, Xj, Aj) = {s/2W sin w, v2~W cosiu, Xj, Aj) the Hamiltonian reads 
H = H + H 1 +H 1 + n^V2Wo (tpCxj. - C'(A 2 )C* ,i) , 



(35) 



where the notations are defined by the formulas ( 24 1 to ( 27 1. This Hamiltonian presents 



an additional term with respect to the previous case, the last term, that represents the 



linear term in Y in Hamiltonian ( 15 1. Then, we define a new set of canonical variables 
as: 



z = w - u = w - (X[ + C(A' 2 )), 
#i =Ai, 
#2 = A 2 , 



Z = W, 

0i = A\ + W, 

6> 2 = A' 2 + C' (X' 2 )W, 



(36) 



where z and $2 are slow angles. After having averaged the Hamiltonian (351 with 
respect to the only fast angle #1, the Hamiltonian reads: 



#1:1 =n v //J0 2 + (c - 



1 + 2V£C(02 




(37) 



+ 71 V 2Zo^/ r jlpcos : 



3.3 



where (2, Z, t?2, O2) represent the averaged variables as in section 

As the term y / /IC'('^2) is always small with respect to the unity, it can be neglected 
in first approximatiorj^] Consequently, we are left with the one degree of freedom 
Hamiltonian: 



K 



(o) 

1:1 



(er — n)Z ■ 



z 2 



+ n V 2Za^/JIp cos z, 



(38) 



which corresponds to the classical problem known as the second fundamental model of 
resonance, or also Andoyer's model ( jHenrard and Lemaitre|1983| ). 

The fixed points of the Hamiltonian system correspond to the slow component of the 
forced orbits (average with respect to the angle fast angle Ai). If we set Q = V 2Z, the 



3 If the Trojan satellite is at L4 or L5, then £'($2) = Vt. 
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coordinates (z, Z) of these fixed points are given by the positive roots of the polynomial 
equation: 



P~ (Q) = Q 3 - 4(ff - n)Q - Anp^pF = if z = 0, 
P + (Q) = Q 3 - 4(cr - n)Q + -inp^TpB = if z = 7T. 



(39) 
(40) 



As P~ (-Q) = —P (Q) , the roots of one of these polynomials are the opposite of the 
roots of the other one. 

By fixing the parameters fj,, p, n and vary the frequency a (in the neighborhood of 
n), there exists a critical value ao such that for a < ctq the average resonant system 
has only one fixed point, which is stable. Then, when a reaches the critical value ctq 



a pitchfork bifurcation with symmetry breaking (see e.g. Lichtenberg and Lieberman 



( [1992| ) occurs and gives rise to two new branches of fixed points. One of these new 
families contains stable fixed points while the other one possesses unstable equilibrium 
points. The bifurcation arises when the discriminant A of the polynomial involved in 
( 39 1 or ( 40 1 is equal to zero, that is 



A — 16 ( n 2 app 2 



16 
27 



(a - n) 3 ) 



= 0, 



165 — 27 pp n 8 — 27 pp n = with 5 = a — n 
Since pp 2 < 4, this equation has only one real root 5q given by: 

1/3 



i - 



pp< 



1/3' 



(41) 
(42) 

(43) 



3n , 



(4 / xp 2 ) 1 / 3 + (M 1/J ) 



l/3s 



In other words, the critical value of a is: 



o-q = n + 5q = ( 1 + 



[(W) 1/3 )»+V /3 ). 



(44) 



Finally, the family of fixed points that exists for all values of a is gi ven by (z\, Z\) = 
(0, Qi/2) where Q\ is the positive root of the polynomial equation (39 1 given by the 
expression 



1 = ( 2np^/pa + 



1/3 



+ I 2npy/pa 



1^) 



1/3 



(45) 



where r\ = 1 if a < o"o and r\ — if a > uq. 

The two others fixed points which coincide when a — cro and split for larger values 
are of the form (z p ,Z p ) = (ir,Qp/2), where the Qp's (p = 2,3), which are the two 
positive roots of P + , can be calculated in the following way. Let us define the complex 
number u by the expression 



u = — 2npy'pa + 



(46) 



As the real part of u is negative and its imaginary part positive, <f>, the argument of u, 
is between ir/2 and 7T. It turns out that the two positive roots of P + read: 



= 4 



= 4 



(j> + 47T 



(47) 
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It is worth mentioning that, as P (— Q) = —P + (Q), the positive root of P can 
also be read: 

4> + 2vr 



a — n 
-4\ / — - — cos - 



for a > (Tq. 



(48) 



3 3 

The expression of the Qj are particularly simple when a — ctq. Indeed, as for this value 
of a the argument <j> is equal to 7r, we have 



Qi = M 



h = Qs = 2\ 



(49) 



Regarding the stability of these three equilibrium points, a straightforward calcu- 
lation shows that the characteristic polynomial of the linearized differential system in 
the neighborhood of (zj, Zj) is given by: 



Vj{h) = h 2 + In^^p^^, 



(50) 



P' being the first derivative of P~ if j = 1 and of — P + if j = 2, 3. 

The polynomial function P~ (Q) being increasing for all values of a, its derivative 
is positive and the roots of Pi are purely imaginary numbers. It turns out that the 
fixed point of coordinates (0, Z\) is linearly stable. The other family is more interesting. 
Indeed, when a = ag, Q2 = Q3 is a double root of P + and consequently its first deriva- 
tive vanishes. This "double" equilibrium point is then degenerated (its eigenvalues are 
both equal to zero). As soon as a becomes larger that gq, the point associated with Qi 
becomes hyperbolic, while the one corresponding to Q3 becomes elliptic. Indeed, the 
polynomial P + having three distinct roots, one negative and the two others Q2 and 
Q 3 positive, the produce (Q 2 ) (Q3) is negative. As Q 3 < Q 2 , ^^-[Ql) > 

anc ^ ^KJ - ^ 3 ) ^ 0> wn i cn proves the assertion stated above. 

A straightforward computation, based on the determination of the generating func- 
tion used to averaging the Hamiltonian (351 with respect to #1, allows us to express 



the fixed points of (38 1 in the original variable (x,X). In particular, we have: 



Q.> 



s'm(zj + Ai) + o I y/Ji 



■O(v^). 



(51) 



In the previous expression, the index j is equal to 1 if a < <jq (only one fixed point) 
and to 1, 2 or 3 otherwise (three fixed points). When a is close to the bifurcation value 



o"0 , according to the expressions ( 44 1 and ( 49 1 , we have 



c 3 (2y/Jj,p)3 sin(zj + Ai) + 0(^3), 



(52) 



with ci = 2 and C2 = C3 = 1. Consequently, the amplitude of the resonant case is 
much larger than in the non-resonant case, where the forced trajectory is of amplitude 



0(y/p) (see formula (21 1 and (33 1). 

Now we illustrate the analytical results by some numerical simulations. Until the 
end of this section, we fixed the free parameters of the model such that : fi = 1 X 10~ 6 , 
n = 1, p = 1, tp = — 7r/3 and v Py o = 0. According to the expression (441 the critical 
value a is then equal to gq — 1.011953. 

Let us first illustrate the evolution of the fixed points of the Hamiltonian (|38|) 



close to this bifurcation value. The quantities Qj , which determine the location of the 
equilibrium points through the formula ( 51 \ are plotted in Fig [2]a for a varying in 
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0.6 




a 



Fig. 2 Bifurcation inside the 1:1 secondary resonance. The upper graph displays the values 
of the Qj (Y-axis) versus a (X-axis). The quantities associated with the stable equilibria are 
represented with solid lines while dotted lines are associated with the unstable equilibrium. 
The frequencies of the libration around the two stable fixed points are plotted in the bottom 
graph. The two vertical dotted lines indicate the values of the parameter a used in the numerical 
simulations presented in Fig. [3] 



the interval (0.96, 1.064) including errj- For a smaller than the critical value the phase 
space has only one equilibrium point denned by Q\ (bold curve). In this domain, Q\ 
increases since the distance of the fixed points from the origin is of order \i 1 and of 
order /i 1 ^ 3 when a is close to oq. After the bifurcation, this equilibrium point moves 
away from the origin. 

For a > o"o, the stable equilibrium associated with Q\ persists while a new pair 
of fixed points appears at a distance of order fi ' from the origin: an unstable point 
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related to Q2 (dotted line) and a stable one associated with Q3 (thin solid line). 
While the unstable point moves away from the origin as a increases, Q3 decreases, 
reaching values of order ^fp, outside the resonance. In other words, the stable fixed 
point associated with Q3 plays the same role after the 1:1 resonance than the one 
connected to Qi before this secondary resonance. The figure [2]b shows the evolution 
of the frequencies /1 and /3 of the two elliptical (stable) fixed points deduced from 
the equation Vj( 



lfj) = 0, where the Vj are defined by the relation (501. When a 



increases, the value of the frequency fi begins to decrease, reaches a minimum for a 
slightly greater than 1, and then increases slowly. The fixed point that is created at 
the bifurcation being degenerated the frequency (dotted curve) starts from zero and 
then increases rapidly along the stable branch. 
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Fig. 3 Phase portraits of the 1:1 secondary resonance in the coordinate system (x,x). The 
figures (a) and (c) show the iterates of the first return map of the unperturbed problem for 
two different values of a chosen in both side of the critical values <jq: u = 1 in figure (a) 
and a = 1.013 in (c). The figures (b) and (d) represent the libration frequency f(x) for the 
corresponding values of a. The two gaps visible in figure (d) correspond to the crossing of 
the separatrices where a logarithmic singularity is encountered. The tangent to the frequency 
curve in this point is consequently horizontal. 



In order to get more details regarding the structure of the phase space in the 
neighborhood of the forced trajectory, we use Poincare sections. To this end, we will not 
consider averaged equations, or expansions in a neighborhood of zero as we have done to 
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carry out the analytical computations, but the complete initial differential equation |8| 
where the temporal evolution of v, r and £ is deduced from |2| to (jSJ. As in the above 
analytical development, we have neglected the terms containing the function £ , the 
function £ has to be constant along the integration. This can be done easily assuming 
that the satellite is at the equilibrium points L4 or L5, that is £ = £ = ±7r/3 and 
£ = f = 0. Consequently, we are left with a one degree of freedom system depending 
periodically on the time (just like in the Keplerian case) and the phase portrait can 
be studied using the first return map defined by the flow of the system at time 2n/n. 
The global phase portrait of this problem is well known, it resembles to the one of the 
7r-periodic pendulum weakly perturbed by a periodic time dependent function (see e.g. 



Wisdom (1987 20041). Here we focus on a neighborhood of x = x = because we are 



interested in a local phenomenon that disturbs the forced orbit when a is close to the 
critical value o"o. 

The Figs. [3^. and[3ji show two different phase portraits of the system in coordinates 
(x, x), restricted to a small area around the origin. The Figj3^i is computed with a — I 



chosen below the critical value ctq ~ 1.011953 (Eq. 44 ), whereas on the bottom plot 



(c) a = 1.013 is chosen slightly above the critical value ag. The phase space (a) shows 
only one elliptical fixed point (a periodic orbit) which is not centered on zero due to 
the non-zero value of Q\. For the phase space (c), a second elliptical point appears at 
the position Q3 because the system has crossed the bifurcation located at ctq. Around 
these points, the trajectories rotate in two different directions and these two dynamical 
states are separated by the separatrix emerging from an hyperbolic point located just 
upper the Q3 fixed point. 

The Figs.[3j3 and[3ji represent the graphs of the frequency map of the system for the 
two selected values of a. More precisely, assuming that the system is regular, the orbit 
of the first return map passing through the point of coordinates (0, x) is quasiperiodic 
and has two fundamental frequencies: n and an independent frequency denoted f(x). 
The frequency map is the map which associated with x the frequency f(x). In prac- 
tice, this frequency is evaluated numerically applying the frequency analysis developed 



by Laskar (1990 19991. If the dynamics of the studied system is regular (integrable 



system), the graph of the frequency map is almost everywhere perfectly smooth. On 
the contrary, its singularities or lack of smoothness indicate the location of chaotic 
zones, or of hyperbolic domains associated with resonances. In addition, when a stable 
island associated with a resonance is crossed, the frequency f(x) remains constant. For 
further explanations and developments about the frequency map analysis, we refer the 
reader to |Laskar| ( |1999[ ) and references therein. 

The frequency /(£) versus x is plotted in the right column of the figure [3] In 
Fig. [3]b, where a = 1, the curve seems perfectly smooth, which indicates the global 
stability of the considered problem (if unstable regions exist, their size is microscopic 
at the scale of the figure). The frequency reaches a local maximum corresponding to 
the fixed point. Its value fits very well to the value predicted by the relation (50 1. The 
two local minima are not directly related to critical points or singular structures of the 
phase space, but by increasing the values of <r, these minima will tend to zero. When 
(7 — (To, one of these minima indicates the location of the fixed point Q2 = Q3 and the 
other one, the location of its separatrix. Notice that, for a = 1, the frequency f(x) is 
always greater that 0.01, that will be of interest when the coorbital perturbation will 
be added to the model. 

For a — 1.013 (Fig. |3]d), the frequency curve is globally smooth, but possesses 
three singularities associated with two gaps (close to x = 0.15 and 0.18) and a cusp (at 
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x — —0.4). These structures correspond to the crossing of the separatrices emerging 
from the unstable equilibrium. Between these singularities, the frequency reaches two 
extrema, the positive one, which coincides with the equilibrium located inside the big 
clockwise libration island, and the negative minimum associated with the fixed point 
lying inside the small counterclockwise libration island. Here again, the agreement 
between the values of the extrema and the values predicted by fx and is better than 
1%. 




-6 -4 -2 2 4 6 8 10 

f/v 

Fig. 4 Chaotic behavior of the libration frequency of the perturbed problem: enlargement of 
the counterclockwise libration island. The X-axis represents the ratio of the frequencies f(x) / v, 
thus the abscissa of the main resonance integer (or half-integer) . The values of x are on the Y- 
axis. This figure highlights the overlap of the lower order resonances and its associated chaotic 
area. The regions where f(x) is constant (vertical segments) are associated with the resonances 

f/u = P en. 



Let us now study the effect of the coorbital perturbation on the rotation inside 
the 1:1 secondary resonance. As the perturbed problem has three degrees of freedom, 
we cannot visualize its dynamics using Poincare sections. But we can still study the 
frequency map of the perturbed problem. Here the situation is very different from the 
non-resonant case studies in section [3li| Indeed, in the non-resonant case, the libration 
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frequency is close to a. It turns out that the ratio of this frequency by the coorbital 
frequency v is of order 1 / ^//Z. Then, the possible resonances between these frequencies 
are of high order and their dynamical influence is negligible. But, as in the resonant case, 
the libration frequency can reach null values, so we can expect to observe interesting 
phenomena that we detail below. 

As the size of the perturbation increases with the amplitude of the variations of £, 
a large amplitude tadpole orbit has been chosen to introduce the perturbation in our 
simulation. Its initial conditions are f (0) = 7r/3, ((0) = 0.0024. Along this orbit, £ varies 
between 24.3° and 162.5° and the coorbital frequency is equal to v = 1.5085 10 -3 . 

In the case of a = 1, the frequency curve of the perturbed problem (not presented 
in the paper) is extremely similar to the one of the unperturbed problem drawn in Fig. 
|3]b. This shows that the addition of the perturbation does not generate any significant 
instability. This is probably due to the fact that the minimum of the libration frequency 
being of about 0.0103, the lower order reachable resonance is given by / = Iv, and the 
system is in some sense protected from the destabilizing influence of the perturbation. 
Consequently, the chaotic zones resulting from these resonances are extremely small. 

On the contrary, the situation is much more interesting when a = 1.013 for at least 
two reasons. First, due to the presence of a hyperbolic fixed point, very thin chaotic 
regions exist near the separatrices of the unperturbed problem. Then we expect that, 
adding one degree of freedom, the size of these chaotic region increases. The second 
reason lies on the fact that, as shown in Fig [3]d, the frequency / approaches zero, 
making possible resonances of low order between / and v. The expected phenomena 
are observed in Fig [4] which presents an enlargement of the graph of the perturbed 
frequency map. The global shape of the frequency curve is similar to the case of the 
unperturbed problem (Fig. [3]d), but new structures are clearly visible. This curve 
presents two very different regions. The first one is composed of three pieces of the 
curve: the part lying above the line x — 0.196, the part below x = 0.025 and the one 
located in the left side of the vertical line f/v= —4. In this first domain, the curve is 
very smooth except for some regions marked by discontinuities generated by resonances 
between f(x) and v. But they are isolated and their associated chaotic areas are not 
observable, even at this scale. 

The second region includes what remains of the two (quasi) horizontal branches 
of Fig. [3]d marking the crossing of the separatrices of the unperturbed problem. The 
regular branches are replaced by vertical segments, whose abscissa is integer, and by 
regions where the points are not well arranged. These structures are mainly generated 
by the presence of low order resonances / = kv, the integer k satisfying —4 < h < 5, and 
of their overlapping. Although these phenomena are clearly visible, they are confined 
to the vicinity of the separatrices of the unperturbed problem. In the parameter space 
that we explored, this should not significantly affect the rotation of the satellite. By 
extension, in the case of the coorbitals, with the parameters listed in Tables [T] and [2] 
the resonances should stay confined in very limited phase space. 



3.4-3 The 2:1 secondary resonance 



According to Section 3.3 the invariant set (y'{t), Y (t), \\ (t), A^(t)) = (0,0, nt + X^ 



(0) 



y/Jint + X'^) corresponds to the forced trajectory which is elliptic (stable) in absence 
of secondary resonances. As its normal frequency is close to a, a doubling period like 
bifurcation arises when a is close to n/2. In order to study this phenomenon, we 
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introduce the canonical coordinates: 

z = 2w - u = 2w - (A' x + C(A 2 )), Z = W/2, 

#1 =Ai, 6>i = A[ + W/2, (53) 

tf 2 =A 2 , 9 2 =A' 2 +('(\' 2 )W/2, 

such that z and $2 are two slow angles. After the substitution of (to, W, \'j,A'j) by the 
new coordinates (z, Z,$j ,Oj) in the Hamiltonian H defined by the relations |24l to 



(26 1, followed by its averaging with respect to #i, the resonant Hamiltonian reads 



K 2 -.i = (2a - n)Z - Z 1 + (3aZ - 4Z 2 ) cosz ^ 

+^n0 2 + ^L{2a-n- 2Z) Z?{& 2 )- 

As in Section |3.4.1| in a first step, we can focus on the first line of the previous 
expression, neglecting the terms depending on (& 2 ,0 2 ). The topology of the phase 
space of this Hamiltonian system whose Hamiltonian is denoted by bifurcates for 
two different values of the parameter a: 

n/„ 3 ^ \ _1 , n A 3 



-0=2 



(l - - v^p) » and a * = \ f 1 + _ vW) ' ^ 55 ^ 



When a < <jq the system associated with K^) does not have any fixed point, except 
the point Z — 0, which corresponds to the singular half-line of the polar coordinates. 
In other words, in cartesian coordinates (V2Z sin z, y2Z cos z), the origin is a stable 
fixed point, and the topology of the phase space is similar to the non resonant case. 
When the first critical value ao is reached, the origin (or the half-line Z = 0) loses its 
stability given rise to a new elliptical fixed points of coordinates (zq, Zq) with 

2(2cr - n) + S^Jipa 
20 = °' Z ° = 4(1 + 2^) ■ (56) 

From the second critical value a^, the central equilibrium point recovers its stability 
while a new unstable fixed point of coordinates , Z^), where 

2(2<r - n) - S^Jipa 
^ = 4(1-2^) ' (57) 

emerges from Z = 0. Its separatrix surrounds the point (zo, Zo) which becomes stable 
again. This kind of bifurcation is described by Murray and Dermott (1999, Chap. 8, Fig. 
8.12) for a different situation corresponding to an internal second-order mean motion 
resonance in the restricted three body problem. One gets the same phase portrait as 
in the book using a coordinate system like (Z cos(z/2), Zsin(z/2)). 

As it is the case for the 1:1 secondary resonance studied in Section[3AT] if the value 
of the parameter a exceeds gq, the frequency can approach zero. If we take into account 
the orbital perturbation, this will give rise to resonances with the frequency u, gener- 
ating wide chaotic regions in the neighborhood of the separatrices of the unperturbed 
problem. But the goal of this paper is not to describe these local phenomena. 



4 In order to simplify the notations, the bars symbols over the variables are omitted in the 
expression of the averaged Hamiltonian. 
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4 Discussion and Conclusion 



In previous sections, we have described the dynamical behavior of the rotational motion 
of the coorbital satellites. Notably, we show that three parameters control the dynamics 
(i) the frequency a related to the dynamical figure of the satellite, (ii) the departure 
from the Lagrangian points that leads to introduce orbital perturbations represented 
by £ and adds a new frequency in the system, and (iii) the combination ^ffip coming 
from the short-period orbital perturbation. 

(i) The coorbital satellites are small satellites of some kilometers (the mean radius 
ranges from 1.3 km for Polydeuces to 89.5 km for Janus). The Cassini spacecraft 
provides high resolution images leading to mapping and cartography these satellites 
with an accuracy of 0.4 km for Polydeuces ( Porco et al.|2007 Thomas 2010 1. The only 



accessible parameters for the figure of the satellites are the long axes ellipsoidal shape. 
From these data, we can compute the moment of inertia by assuming that the satellites 
are homogeneous and ellipsoidal 

(58) 
(59) 
(60) 

where A, B, C are the principal moment of inertia and a, 6, c are the fitted long axis 



A = 


f" 2 


+ c 2 : 


B = 


M, 2 


+ a , 


C = 


M 2 


+ b 2 



ellipsoidal shape coming from shape models (Thomas 20101. By looking the recent 
images of Helene (Cassini image N00172886.jpg, NASA/JPL/Space Science Institute), 
it is clear that it is a rough hypothesis. We use the data provided by |Thomas| |2010 ) 
that compiled the recent shape models of these satellites and we compute the triaxiality 
(B — A)/C that is equal to 



B - A a 2 -b 2 



C 



(61) 



so the main uncertainties come from the equatorial axes of the shape. The resulting 
values of a = n^3(B - A)/C are shown in Table ^ To give an idea of the error 
due to the ellipsoidal approximation, it is interesting to review the determination of 
the shape of Janus and Epimetheus. Tiscareno et al.| ( |2009[ ) provided the moments of 



inertia of both satellites from shape model without the assumption of the ellipsoidal 
shape and deduced a equal to 8.51, and 4.94 rad/day instead of 5.62, and 3.85 rad/day 
for the ellipsoidal shape model. Clearly, the ellipsoidal shape model gives an error of 
the order of 33 %. The Table [2] provides the different values of a for the coorbitals by 
assuming an ellipsoidal shape. We note that Telesto, Calypso, and Polydeuces present 
a departure from the 1:1 secondary resonance of 3.3, 12.0, 23.2 %, whereas Helene 
shows 19.0% departure from the 1:2 secondary resonance. Consequently, by assuming 
an uncertainty of 30%, these satellites could be in or close to secondary resonance 
and more data are necessary to determine their figure. Basing on the values given in 
the tabels [T] and [2] we estimate the amplitude of the forced orbit for each of the six 
coorbital satellites. To this end, we assume that the rotation of the satellites is not 
in secondary resonance, and use the expression ( |21[ | to derive the amplitude of the 
libration with respect to the line Saturn-satellite Ag = 2p^fjln 2 (n 2 — o 2 )~ l , which 
is given in the fifth column of the table [2] The amplitude of the forced orbit in an 
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unif ormly rotating reference frame (Tab. |2| last column) derived from the expression 



(22 l is given by A r = 2p^/Jlo (n 



(ii) The departure of the trajectories from the Lagrange points leads to long-period 
oscillations in the orbital motion characterized by (. The perturbation of the orbital 
motion influences the rotational motion of the satellites by varying the gravitational 
torque coming from Saturn expressed in jSJ, through r and v. This is an indirect effect 



acting on the satellite. It has been shown in the equation (21 1 that the amplitude of 
the forced trajectory is almost the same than in the Keplerian case but its phase differs 
by the ^-periodic function (,(yi). For the a values listed in the Table [5J the influence of 
£ is small for Calypso, Telesto, and Helene. However, for large tadpole orbits such as 
Polydeuces and for horseshoe orbits (Janus and Epimetheus) the influence is large. In 
addition, for a leading to secondary spin-orbit resonance, the presence of the additional 
frequency v causes new resonances and generates chaotic areas. 

(iii) The orbital motion of the satellites is described through the Erdi' s formal 



ism (Erdi 1977) and the initial values come from the Horizons ephemerides ( Giorgini 
|et al.|[l996[ ). The Horizons ephemerides provide a recent orbit of the satellites based 
on additional images from Cassini spacecraft ( Jacobson et al.|2008 l. In our model the 



parameter y/Jip controls the amplitude of the orbital short period variations as the 
eccentricity does for Keplerian motions. Consequently, p is directly adjusted to the 
mean eccentricity of the satellite. The eccentricity of the coorbitals from the Horizons 
ephemerides varies with time on short and long periods due to several effects. First, the 
eccentricity of the principal satellites is not equal to zero. The eccentricities of Dione 
and Tethys are 0.0022 and 0.0001, respectively. Erdi (1977) shows that in this case, the 
Equation (2) will present an additional periodic term at the orbital frequency leading 
to an additional term of frequency n in the solution. Moreover, due to the interactions 
between the satellites, the longitude v will present secular terms, long periodic compo- 
nent, as illustrated for Janus and Epimetheus (Robutel et al. 20111. All these effects 



could be easily introduced in the present formalism by adding periodic terms in the 
equation of motion. On contrast, the interactions with Enceladus and Mimas that are 
in 2:1 resonance with Dione and Tethys, and therefore with the coorbitals, are more 
complicated to model and to describe. 



Satellite 


n 


a 


y/Up 


A S 






(rad /day) 


(rad/day) 


(dcg.) 


(deg.) 


Polydeuces 


2.3 


1.86 


0.0196 


6.5 


4.2 


Helcnc 


2.3 


1.42 


0.0078 


1.5 


0.6 


Telesto 


3.3 


3.22 


0.0006 


1.4 


1.4 


Calypso 


3.3 


2.97 


0.0008 


0.5 


0.4 


Janus( a ) 


9.03 


4.94 


0.0068 


1.1 


0.3 


Epimetheus^ 


9.03 


8.51 


0.0098 


10 


8.9 



Table 2 Amplitude of the forced orbit. The mean motion n and the parameter <r (libration 
frequency at the center of libration) are given in the second and third columns. The fourth 
column corresponds to the values of the parameter sjpp that governs the variations of the 
orbital motion at the frequency n. A$ is the amplitude of the forced orbit measured with 
respect to the line Saturn-satellite, while A r represent the amplitude of the short-period term 
of the forced tra jectory in uniformly ro tating frame. The symbol (a) indicates that the value 
of a comes from Tiscarcno ct al. ( 2009| . 
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In conclusion, we have studied the impact of the coorbital resonance on the ro- 
tational motion of Telesto, Calypso, Helene and Polydeuces, in the planar problem 
(orbital and equatorial inclination are neglected). We have shown that the main pa- 
rameters acting on the rotational motion are (a, £, y/pp) where a is related to the 
dynamical figure of the moons, £ is the secular perturbation of the orbit and y/pp is 
the departure from the Lagrangian points. The uncertainties in the dynamical figure 
of the moons range the possibility that they are close or in secondary resonances. Such 
secondary resonances, coupled with the secular motion of the moons due to the coor- 
bital resonance, increase the possibility to reach supplementary resonances and these 
satellites are inclined to exhibit a rich dynamics. The determination of the dynamical 
figure of these satellites by space observations will bring precious information for the 
dynamics of these moons. 



Appendix 



(19771 and the equations \2\ to |5 



The link with the expressions given by Erdi 
not straightforward and some cl ues are given in the following lines. The x tn formula 



Erdi 



( 1977 1 and by (x) if it comes from the present 



is denoted by [x] if it appears in 
paper. 

From the relations [5a] and [19] truncated at first order in e = ^/p, we deduce that 
r = 1 + (pi cos(<p + tpi) — §/i), while [5b] and [20] give 8 — (p — 2^/ppi sin(<^ + ?/'i) + 
y/Jlqi (the polar angle 9 is denoted by v in our paper). The functions tp, /i, pi and V>l 
depend a priori on v which is a slow time defined by v — */Ji[v — «o), v being the true 
anomaly of the massless body (v = nt in the Circular Restricted Three Body Problem 
(C.R.T.B.P.)). As from [27a] and [27b], the functions pi and ip\ are constant (ei = 
in the C.R.T.B.P.), it only remains to give an explicit expression for tp and /i . From 
[7a] we have tp — v + ^fp J f\(v)dv and the relation [23] links f\, and consequently tp, 
to £i (£ in the present paper) that satisfies the differential equation [24] (equation (5) 
in this text). By some changing of notations, we get the approximations ([2| to |5|. 

Here, we have extracted from Erdi's paper only the necessary relations to get a 
simple orbital model of a Trojan satellite. However, Erdi's theory goes far further since 
it gives asymptotical expansions to second order in ^fp within the framework of the 

1 1977) summarizes 



planar elliptic restricted three body probl em. The sec tion 7 of 
these results. A spatial theory is given in Erdi (19781. 



Erdi 
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